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ABSTRACT 

We present optical photometry of 16 transits of the super-Earth GJ 1214b, allowing us to refine 
the system parameters and search for additional planets via transit timing. Starspot-crossing events 
are detected in two light curves, and the star is found to be variable by a few percent. Hence, in 
our analysis, special attention is given to systematic errors that result from star spots. The planet- 
to-star radius ratio is 0.11610 ± 0.00048, subject to a possible upward bias by a few percent due to 
the unknown spot coverage. Even assuming this bias to be negligible, the mean density of planet 
can be either 3.03 ± 0.50 g cm~^ or 1.89 ± 0.33 g cm"-^, depending on whether the stellar radius is 
estimated from evolutionary models, or from an empirical mass-luminosity relation combined with 
the light curve parameters. One possible resolution is that the orbit is eccentric (e ~ 0.14), which 
would favor the higher density, and hence a much thinner atmosphere for the planet. The transit 
times were found to be periodic within about 15 s, ruling out the existence of any other super-Earths 
with periods within a factor-of-two of the known planet. 

Subject headings: stars: planetary systems — planets and satellites: individual (GJ 1214b) — stars: 
individual (GJ 1214) — techniques: photometric 



1. INTRODUCTION 

The recently discovered planet GJ 1214b (Cliarbon- 
neau et al. 2009) is the smallest known cxoplanct for 
which the mass, radius, and atmospheric properties arc 
all possible to study with current technology. It is there- 
fore a keystone object in the theory of planetary inte- 
riors and atmospheres, and has been welcomed as the 
harbinger of the "era of super-Earths" (Rogers & Sea- 
gcr 2010). 

The planet's discoverers estimated the mass and radius 
of GJ 1214b to be6.55±0.98Me and 2.68±0.13 i?®, giv- 
ing a mean density of 1.87±0.40 g cm~'^ (Charbonneau et 
al. 2009). This is such a low density that it would seem 
impossible for the planet to be solid with only a thin, 
terrestrial-like atmosphere. Rather, it seems necessary 
to have a thick gaseous atmosphere, probably composed 
of hydrogen and helium but possibly also of carbon diox- 
ide or water (Charbonneau et al. 2009, Rogers & Seager 
2010, Miller-Ricci & Fortney 2010). 

In this paper, we report on observations of additional 
transits of GJ 1214b. As in other papers in the Transit 
Light Curve (TLC) series (Holman et al. 2006, Winn et 
al. 2007), one of our goals was to refine the basic sys- 
tem parameters, and thereby allow for more powerful 
discrimination among models of the planet's interior and 
atmosphere. Another goal was to check for any non- 
periodicity in the transit times, as a means of discover- 
ing other planets in the system, through the method of 
Holman & Murray (2005) and Agol ct al. (2005). Super- 
earths have frequently been found in pairs or even triples 
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in compact arrangements (Lo Curto ct al. 2010), and it 
would be interesting to know if GJ 1214b is another such 
example. 

This paper is organized as follows. Section [2] describes 
the observations and data reduction. Section [3] presents 
the light curve model, taking into account the effects of 
starspots. Section|3]discusses the method by which we es- 
timated the model parameters and their confidence inter- 
vals. Section [5] discusses the results for the planet-to-star 
radius ratio. Section|6]presents two different methods for 
determining the stellar radius (and hence the planetary 
radius), which give discrepant results. Some possible 
resolutions of this discrepancy are discussed. Section [7] 
presents our analysis of the measured transit times, and 
constraints on the properties of a hypothetical second 
planet. Finally, in Section [51 we discuss the implications 
of our analysis on the understanding of GJ 1214b and 
more broadly on M dwarf transit hosts. 

2. OBSERVATIONS AND DATA REDUCTION 

Our data were gathered during the 2009 and 2010 ob- 
serving seasons. Thirteen transits were observed with 
the 1.2m telescope at the Fred L. Whipple Observatory 
(FLWO) on Mount Hopkins, Arizona, using Keplercam 
and a Sloan z' filter. The blue end of the bandpass was 
set by the filter (approx. 0.85 /im) and the red end by 
the quantum efficiency of the CCD (ranging from nearly 
100% at 0.75 Aim to 10% at 1.05 /im). The first two of 
the FLWO transits were already presented by Charbon- 
neau et al. (2009); those data have been reanalyzed for 
this work. Another three transits were observed with the 
6.5m Magellan (Baade) telescope at the Las Campanas 
Observatory in Chile, with the MagIC and IMACS cam- 
eras and a Sloan r' filter (approx. 0.55-0.70 /.im) 

We used standard procedures for bias subtraction and 
flat-field division. Circular aperture photometry was per- 
formed for GJ 1214 and several comparison stars (5-6 for 
the FLWO images, and 2 for the Magellan images). Cir- 
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cular annuli centered on each star were used to estimate 
the background level. The flux of GJ 1214 was divided 
by the sum of the fluxes of comparison stars, to produce 
a relative flux time series for GJ 1214. The radii of the 
apertures and annuli were varied in order to seek those 
values that minimize the scatter in the residuals, while 
also producing a minimal level of time-correlated noise, 
as determined by the wavelet-based algorithm described 
by Garter & Winn (2009). Because the amplitude of the 
correlated component was always < 1% of that of the 
white component, in what follows we assume the noise 
to be uncorrelated in time. Each measurement was asso- 
ciated with a time stamp given by the midexposure time 
expressed in the BJDtdb system (Eastman et al. 2010). 

Fig. [T] shows the light curves, after correcting for dif- 
ferential airmass extinction as described in § 111 Fig. [2] 
gives a clearer view of two of the Magellan light curves 
(epochs 5 and 260), which show evidence for anomalous 
brightening events during the transit. The data points 
corresponding to these events are shown with unfilled cir- 
cles, rather than gray circles. Correlations were sought 
between the flux observed during these events, and var- 
ious image statistics such as the median X and Y pixel 
position, the shape parameters of the point spread func- 
tion, and the background level, but none were found. 
We interpret each of those events as the passage of the 
planet in front of a dark starspot, as has been observed 
by Rabus ct al. (2009) and others for different transiting 
systems. 

3. TRANSIT LIGHT CURVE MODEL 

The evidence for starspots on GJ 1214 caused us to 
make a few changes to the usual models for transit pho- 
tometry. Gool starspots on the visible hemisphere of the 
star can produce two different effects: when they are on 
the transit chord, they cause brightening events in the 
light curve such as those we observed; and when they 
are away from the transit chord, they reduce the over- 
all brightness of the star and increase the fractional loss 
of light due to the planet. The latter effect causes the 
planet to appear larger than it is in reality, by producing 
a deeper transit. 

3.1. Loss of light due to starspots 

During a transit, the flux received from a spotted star 
can be written 

F{t)=F,[l-e{t)]-AF{t), (1) 

where i<o is the flux that would be received from an 
unspotted and untransited star, e{t) is the fractional loss 
of light due to starspots, and ^F{t) is the flux that is 
blocked by the planet. The factor e{t) changes on the 
timcscalc of the stellar rotation period (^ months), while 
Ai^(i) changes on the much shorter timescale of the tran- 
sit (~ hours). 

It is convenient to normalize the light curve by the flux 
measured just outside of the transit, 

where e is taken here to be a constant for the duration 
of the transit. The quantity AF{t)/FQ is the fractional 
loss of light due to the transit, which can be calculated 



in the usual way. In this work we use the formulas of 
Mandel & Agol (2002) for a quadratic limb- darkening 
law. Thus, our model adds one new parameter e spcciflc 
to each transit. A similar model was used for Gorot-2b 
by Gzcsla et al. (2009). 

With infinite photometric precision, e could be deter- 
mined from a single transit light curve based on the slight 
differences in the ingress or egress data as compared to 
a spot-free model. However, the differences are typically 
no greater than a few parts per million, and will be diffi- 
cult to detect in practice (see the right panel of Fig. [3]) . 
Consequently, there is a strong degeneracy between e and 
the planet-to-star radius ratio Rp/Ri,, in the sense that 
increasing either parameter results in a deeper transit 
with no other observable changes. If the spots produce 
an overall darkening (e > 0) then the measured transit 
depth S (corrected to remove the effects of limb darken- 
ing) can only be used to place an upper bound on the 
radius ratio: 



^ = 7^(137) <V-S. (3) 

If the spot coverage is time-variable, then the transit 
depth will vary from event to event. One may hope that 
the shallowest observed transit corresponds to a time 
when the visible disk was nearly free of spots, thereby 
allowing the planet-to-star radius ratio to be measured 
with minimal bias. However, there is no guarantee that 
the spot coverage ever falls to zero. 

3.2. Spot crossings during transit 

When the planet transits a cool spot, the fractional 
loss of light temporarily decreases, and a brightening is 
observed in the light curve. The duration and ampli- 
tude of the brightening event depends on the size of the 
spot and the intensity contrast with the unspotted stel- 
lar photosphere. If the spot is approximately the same 
size as the planet, and is completely occulted, then the 
brightening event will have a duration of approximately 
twice the transit ingress duration. For a circular spot, 
the amplitude of the brightening event is 

where Rs and Is arc the radius and mean intensity of the 
spot, and Z^, is the mean intensity that the photosphere 
would have at that location in the absence of spots. 

Rather than modeling the spots (e.g., Wolter et al. 
2009) and fitting for the brightening events, we chose the 
simple approach of identifying the brightening events vi- 
sually, and assigning those data zero weight in the fitting 
process. For the Magellan data, the identification could 
be performed without much ambiguity, and the "excised" 
data are shown in Figs. 1 and 2. For the FLWO data, no 
brightening events were obvious to the eye. This is prob- 
ably due to a combination of the lower precision of the 
FLWO data, and the weaker contrast between the spots 
and the photosphere at the longer wavelengths of the 
FLWO observations (z' band as compared to r' band). 

In particular, the anomalies we observed at epochs 
5 and 260 have a duration of approximately twice the 
ingress duration, and therefore the typical spots are likely 
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Fig. 1. — GJ 1214 transit light curve photometry. Each panel shows transit light data (after correcting for airmass variation), the best- 
fitting model, and the residuals (observed — calculated, in parts per thousand). Unfilled circles represent the suspected "spot anomalies" 
that were assigned zero weight in the fitting process. The data from epochs 5, 212, and 260 are based on r'-band observations with the 
Magellan Clay 6.5m telescope. The other data are based on 2:'-band observations with the FLWO 1.2m telescope. 



to be comparable in size to the planet. The observed am- 
plitudes of 2-3 ppt (parts per thousand) give a constraint 
on the spot size and the intensity contrast. The lowest 
possible intensity contrast is obtained for the largest pos- 
sible spots, with Rs ~ Rp. In that case, assuming the 
unspotted photosphere and the spots to be described 
by blackbody spectra, the spot temperature would be 
Tg « 2900 K as compared to the photospheric temper- 
ature of 3026 K (Charbonneau et al. 2009). Had we 
observed those same events in the z' band, the brighten- 
ing amplitude would have been wl.S ppt, just below the 
la level of the noise in each data point. Smaller spots 
would need to have a greater intensity contrast, leading 
to a greater temperature difference and an even smaller 
predicted amplitude in the z' band. 

In this light it is not surprising that we did not detect 
similar starspot events in our z'-band data. In what fol- 
lows we assume that the z'-band data are unaffected by 
spot-crossing events. Of course it is possible that smaller 
or less luminous starspots were transited during our ob- 
servations, and that the transit depths are "filled in" to 
some degree by brightening events that cannot be identi- 
fied visually in the light curves. However, the reasoning 
in the preceding paragraph suggests that this effect is 
minor. Furthermore, we did not find any evidence for 
time correlated noise in our model residuals. 



3.3. Using transit durations to estimate the radius ratio 

As discussed, the effect of dark spots on the visible disk 
of the star is to increase the fractional loss of light during 
the transit, causing the planet to appear larger than it is 
in reality. This may be counteracted to some degree by 
the "filling in" of the transit due to spot occultations, if 
those events cannot be recognized and excised from the 
data. This caused us to wonder what one could still learn 
about the planet if the transit depth cannot be trusted. 

The transit light curve provides three primary observ- 
ables: the total transit duration Ttot, the ingress or egress 
duration r, and the depth 6. If one also knows the orbital 
period P, eccentricity e, and argument of pericenter w, 
then one can translate the transit observables into three 
parameters more directly related to the star-planet sys- 
tem. One possible choice of three "physical" parameters 
is the planet-to-star radius ratio Rp/R^,, the orbital in- 
clination i, and the mean stellar density Pi, (see, e.g., 
Carter et al. 2008 or Winn 2010). Usually, the value of 
Pi, derived in this way is used to improve the character- 
ization of the host star (see, e.g., Sozzetti et al. 2007, 
Holman et al. 2007). However, if one is willing to adopt 
a value for based on external information (such as the 
star's observed luminosity and spectrum), then the two 
observables Ttot and r are sufficient to specify the other 
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Fig. 2. — Close-ups of the two light curves with suspected spot anomalies. Both light curves arc based on r'-band observations with the 
Magellan Clay 6.5m telescope. The plotting conventions are the same as in Figure FlT 



two physical parameters i and Rp/R^,. 

Hence one may derive Rp/Ri, based on Ttot and r 
alone, without any reference to 6. The result of straight- 
forward algebra is 



total 
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where G is the gravitational constant. For cases where 
star spots are likely to be present, this "duration-based" 
method of deriving i?p/i?* is a useful alternative to the 
usual "depth-based" method, and it may be more ac- 
curate in situations when pi, is well constrained. The 
duration-based method is immune to the bias that is 
caused by starspots outside of the transit chord, although 
one must still be wary of starspot anomalies that oc- 
cur during ingress or egress. As we will describe, for 
GJ 1214 we considered both depth-based and duration- 
based methods to estimate i?p/i?*. 

4. TRANSIT LIGHT CURVE ANALYSIS 

We modeled the data that are plotted in Fig. [T] using 
the approach represented by Equation ([TJ . Each transit 
was assigned a parameter e, representing the fractional 
loss of light due to spots on that night. Since a degener- 
acy prevents the determination of all of the e values along 
with Rp/Ri,, we fixed e = for the transit with the shal- 
lowest depth (epoch 236, UT 2010 June 5) and allowed 
all the others to be free parameters. This is equivalent 
to assuming that the visible stellar disk had no spots on 
that particular night. The results can be scaled appro- 
priately for any assumed value of the spot coverage on 
that night. 

We used the Mandel & Agol (2002) equations to cal- 
culate AF/ Fq, assuming a quadratic limb-darkening law. 
The orbit was assumed to be circular. The time of con- 
junction tc for each transit was allowed to be a free pa- 
rameter, but for the limited purpose of computing the 
scaled orbital distance a/i?*, we assumed the orbital pe- 



riod to be 1.5803925 days (Charbonneau et al. 2009). We 
also allowed the out-of-transit magnitude on each night 
to be a linear function of airmass, in order to allow for 
color-dependent differential extinction between GJ 1214 
and the bluer comparison stars. 

All together, there were 70 model parameters: Rp/Ri,, 
a/Ri,, i, Ur and Vr (the r'-band limb darkening coef- 
ficients), Uz and Vz (the z'-band limb darkening coef- 
ficients), Ej (for j = 1 to 15), tcj (for j = 1 to 16), 
and two constants specifying the linear airmass correc- 
tion for each of the 16 transits. We determined the pos- 
terior probability distributions for these parameters with 
a Markov Chain Monte Carlo (MCMC) algorithm, using 
Gibbs sampling and Metropolis-Hastings stepping (see, 
e.g., Tegmark et al. 2004, Holman et al. 2006). The like- 
lihood was taken to be exp(— x^/2), with 



[Fk-ntk)] 



(6) 



where F is the measured flux, is the calculated flux, 
and tTfe is root-mean-square of the residuals from the best- 
fitting model specific to each night. The data from the 
suspected spot-induced anomalies during epochs 5 and 
260 were excluded from this sum. 

Figures [T] and [2] show the minimum-x^ solution (black 
curve) and the residuals for each night. Figure S] shows 
the time-binned composite z' and r' light curves, with 
2 minute sampling. The data were combined after cor- 
recting the depths for spot variability [such that F i-> 
1 — (1 — F){1 — Cj)]. Table [1] gives the results for each 
parameter, based on the 15.8%, 50%, and 84.2% values 
of the cumulative posterior distribution for each param- 
eter, after marginalizing over the other parameters. This 
table also gives some other parameters based on subse- 
quent steps in the analysis, described in the sections to 
follow. 

In order to assess any possible transit duration varia- 
tions, we also performed a second analysis, in which the 
orbital inclination i was allowed to be a free parameter 
specific to each night. Table [2] gives the transit depth 
[{Rp/Ri,)^ /{I — ej)], time of conjunction, and total tran- 
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Fig. 3. — The effect of untransited starspots. Left. — The Hght curves of a planet transiting a star with (dashed hne) and without (solid 
line) a starspot on the visible disk [e = 5% in Eqn. ((ij]. In both cases the out-of-transit flux was set equal to unity. Right. — Residuals 
between the light curves and the best-fitting spot-free light curve. 
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Fig. 4. — GJ 1214 composite transit light curves, in the ^'-band (top) and r'-band (bottom), after merging all the data and resampling 
into two-minute bins. The black curves show the best-fitting model light curves. The rms residual is 270 ppm for the r'-band light curve, 
and 560 ppm for the 2:'-band light curve. 

transit depth, which is beneath our typical measurement 
precision of 5%, and hence it is not surprising that no 
such trend was observed. 

Based upon this finding, one way to estimate the 
planet-to-star radius ratio is to assume that spots have 
not significantly affected the transit depth. We may then 
calculate the radius ratio based on the variance- weighted 



sit duration for each epoch. The depths are plotted as 
a function of epoch in Figure [S] Figure [5] shows the du- 
rations, and the residuals from a linear ephemeris, as a 
function of epoch. 

5. THE PLANET-TO-STAR RADIUS RATIO 

Despite our concerns about the effects of starspots on 
the light curves, we did not find any significant vari- 
ation of the transit depth or and duration with time, 
and we did not find any significant departure from a 
linear ephemeris. In addition, the out-of-transit flux of 
GJ 1214 measured directly from the images was constant 
to within a few percent over the course of our observa- 
tions (see Figure [7]). There was apparently a decline by 
1-2% between 2009 and 2010, followed by a gradual 1-2% 
rise throughout our 2010 observations. This should have 
been accompanied by 1-2% variations in the measured 



average of the 16 transit depths 5j = (i?p/i?^)^(l 



16 



Rp 1 V 3j 

= 0.11610± 0.00048. 



(7) 



where the average is taken over all epochs j. This result 
for the radius ratio is slightly smaller than that reported 
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TABLE 1 
System Parameters of GJ 1214 






Parameter 


Value 




Comment 


Transit Parameters 


Orbital period (day) 

Transit ephemeris (BJDtdb) 

Planet-to-star radius ratio, Rp/Ri, 

Orbital inclination, i (deg) 

Scaled semimajor axis, a/R^, 

Total transit duration, T-potal (min) 

Transit ingress duration, r (min) 

Starspot-unaffected quantity, (p,/pQ)~'^''^(_Rp/i?*) 

r' quadratic limb-darkening coefBcients, {ur,Vr) 

z' quadratic limb-darkening coefficients, (u^^Vz) 


1.58040482 ± 0.00000024 
2,454.980.7487955 ±0.000045 
0.1161 ±0.0048 

90+00 

14 71+0.37 

rr, 70-1-0.49 

0^- '•5-0.35 

5.80_o;45 

0.017 ±0.001 

(0.49 ±0.07, 0.48 ±0.22) 

(0.17 ±0.12, 0.64 ±0.30) 




a 
a 
b 
c 
c 

d 


Orbital Parameters 


Radial velocity semiamplitude, (m/s) 
Orbital eccentricity, e 

Stellar projected rotational velocity, v&vai (km/s) 


12.2 ± 1.6 


< 2.0 




c 

fixed 
e 


Stell 


ar Photometric Parameters & Parallax 






J 
H 
K 

Parallax, vr (mas) 


9.750 ± 0.024 
9.094 ± 0.024 
8.872 ±0.020 
77.2 ±5.4 




e 
e 
e 
e 


Stellar ilfe Planetary Parameters 


Stellar mass, M* {Mq) 

Stellar radius, Ri, {Rq) 

Stellar mean density, p* (g cm~^) 

Stellar effective temperature, T^g (K) 

Stellar surface gravity, log (cgs) 

Stellar metallicity ([Fe/H]) 


Method A| Comment 


Method Bt 


Comment 


0.157 ±0.012 
0.210 ±0.007 

24.1 ±1.7 f 
4.99 ± 0.04 


0.156 ±0.006 
0.179 ±0.006 
38.4 ±2.1 
3170 ± 23 
5.12 ± 0.01 



fixed 



Planetary mass, Mp (Mgj) 
Planetary radius, Rp (-R®) 
Planetary mean density, pp (g cm~^) 
Planetary surface gravity, Qp (m s~-^) 



6.45 ±0.91 
2.65 ±0.09 
1.89 ±0.33 
9.0 ± 1.5 



6.43 ± 0.86 
2.27 ±0.08 
3.03 ± 0.50 
12.2 ± 1.9 



Note. — fParameters calculated assuming an empirical lum inosity-mass relationship (see ^ l6.1l for details). fParameters constrained to 
the Baraffe et al. (1998) stellar evolution isochrones (sec 16.21 for details), a = Determined from a linear fit to the midtransit times listed 
in Table [2] b = Determined by taking the square root of the variance weighted average of the depths listed in Table [2] (see § \5\ for details) . 
c = Assuming epoch 236 coincides with the spot-free stellar surface, d = As calculated with Eqn. [5] c = Charbonneau et al. (2009). f 
= Estimated from transit parameters as p* = (37r)/(GP^)(a/_R,)'^. g = Estimated from transit parameters as (2n/ P)Kir{a/ Rp)^ {sini)~^ 
(Southworth et al. 2008). 



TABLE 2 

Epoch specific transit light curve results 



Epoch Depth, (iip//?j,)2(l - e)"i Midtransit time (BJDtdb) Duration (min) 






0, 


.01332 ±0, 


,00057 


2, 454, 980.74857 ± 0.00015 


52, 


.30 ±0.56 


2 


0, 


.01388 ±0, 


.00059 


2, 454, 983.90982 ± 0.00016 


52, 


.70 ± 0.48 


5 


0, 


.01355 ±0, 


.00059 


2, 454, 988.650808 ± 0.000049 


52, 


.84 ± 0.29 


14 


0, 


.01333 ±0, 


.00063 


2, 455, 002.87467 ± 0.00019 


52, 


.59 ± 0.42 


183 


0, 


.01384 ±0, 


.00059 


2, 455, 269.96299 ± 0.00016 


52, 


.37 ± 0.46 


195 


0, 


.01281 ±0, 


.00081 


2, 455, 288.9282 ±0.0011 


52, 


.64 ± 6.23 


200 


0, 


.01345 ± 0, 


.00067 


2, 455, 296.83013 ± 0.00023 


52, 


.40 ± 0.73 


212=^ 


0, 


.01416 ±0, 


.00066 


2,455,315.79485 ±0.00023 


52, 


.76 ± 0.65 


212'> 


0, 


.01380 ±0, 


.00061 


2, 455, 315.794693 ± 0.000080 


52, 


.82 ±0.31 


214 


0, 


.01303 ±0, 


.00062 


2, 455, 318.95523 ±0.00017 


51, 


.61 ±0.53 


236 


0, 


.01271 ±0, 


.00055 


2, 455, 353.72387 ±0.00018 


52, 


.44 ± 0.66 


238 


0, 


.01377 ±0, 


.00058 


2, 455, 356.88495 ± 0.00015 


52, 


.59 ± 0.42 


243 


0, 


.01313 ±0, 


.00059 


2, 455, 364.78700 ± 0.00015 


52, 


.69 ± 0.33 


250 


0, 


.01299 ±0, 


.00053 


2, 455, 375.84997 ± 0.00013 


52, 


.76 ±0.31 


255 


0, 


.01376 ±0, 


.00060 


2,455,383.75205 ±0.00013 


52, 


.56 ±0.31 


260 


0, 


,01414 ±0, 


.00061 


2, 455, 391.654105 ± 0.000059 


52, 


.72 ± 0.36 



Mean= 0.01348 ± 0.00011 

^ Sloan z' 
Sloan r' 

Variance weighted mean. 
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Fig. 5. — Transit light curve depth as a function of observation epoch (using values tabulated in Table[2]l. 
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Fig. 6. — Transit durations and timing residuals. Left. — Total transit duration versus observation epoch. Right. — Differences between 
the measured times of conjunction, and the best-fitting model assuming a constant period. Only the complete transits presented in this 
work (solid circles) were used to derive the model parameters. Partial transits (open circles) were excluded from the fit. The other plotted 
data (asterisks) are from Charbonneau et al. (2009). 
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Fig. 7. — Out-of-transit variability of GJ 1214. Plotted is the flux (measured out-of-transit) of GJ 1214, relative to four different 
comparison stars. 



in the discovery paper (0.11620 ± 0.00067; Charbonneau 
et al. 2009), but the two results are consistent to within 
the uncertainties. 

Due to time variability of starspots at the 1-2% level, 
this estimate is subject to a bias of a few percent, in the 
sense that the planet-to-star radius ratio may actually be 
a few percent smaller than 0.11610. An even more con- 
servative stance would recognize that wc cannot exclude 
even larger effects due to the time-independent compo- 
nent of the starspot coverage (if, for example, the poles 
of the star were persistently darker than the rest of the 
photosphere). In that sense we can only place an upper 
bound on the radius ratio: i?p/i?* < 0.1161 at 95% con- 
fidence (based on the marginalized posterior distribution 
produced by the MCMC algorithm). 

It should also be noted that the r'-band transit depth 
(1.383% ± 0.035%) was found to be shghtly larger than 
the z'-band transit depth (1.340% ± 0.017%), although 
the difference is only at the l-2a level. A deeper r'- 
band transit is expected if cool starspots are affecting 
the results. In contrast, models of the atmosphere of 
GJ 1214b generally predict that the transit depth should 
have the opposite wavelength dependence, with a deeper 
z'-band transit (see, e.g.. Fig. 3 of Miller- Ricci & Fort- 
ney 2010). Those who would attempt to attribute any 
slight wavelength-dependence of the transit depth to the 
planetary atmosphere should beware of the confounding 
effects of starspots. 

6. THE RADIUS OF THE PLANET 

In order to understand the structure and atmosphere 
of GJ 1214b, we want to know its radius, rather than 
just the planet-to-star radius ratio. This requires some 
externally derived estimate of the stellar radius, or mass, 
or both. As we will describe, we have investigated two 



different pathways to the planetary radius, relying on dif- 
ferent assumptions, and found them to give inconsistent 
results. The results of both of these methods, including a 
number of derived stellar and planetary parameters, are 
given in Table 2. This inconsistency had already been 
noted by Charbonneau et al. (2009), but here we delve 
further into the details and discuss possible resolutions. 

6.1. Method A: Empirical mass-luminosity relation, and 
transit-derived stellar mean density 

In the first method, the stellar mass is estimated based 
on its observed luminosity (parallax, and apparent mag- 
nitude). Then, the stellar radius is found by combining 
the stellar mass with the mean stellar density p^, derived 
from the transit light curve. 



3tt 

GP^ \R* 



1 



Mp 



(8) 



which, for Mp ^ M^,, is only a function of the 
photometrically-determined parameters a/Ri, and P. 

We begin with the ii'-band mass-luminosity function 
of Delfosse et al. (2000). For GJ 1214, with a par- 
allax TT = 77.2 ± 5.4 mas and apparent K magnitude 
8.782 ± 0.020, we used the polynomial function given by 
Delfosse et al. (2000) to estimate = 0.157±0.012 Mq, 
where the uncertainty is based only on the propagation 
of errors in tt and K. There should be additional uncer- 
tainty due to the intrinsic scatter in the mass-luminosity 
relation, but this intrinsic uncertainty was not quanti- 
fied. Charbonneau used the same method and reported 
= 0.157 ±0.019 Mq. 

Next, using a/Ri, = 14.71 ± 0.35, as derived from our 
ensemble of light curves, we found the mean stellar den- 
sity to be = 24.1 ± 1.7 g/cia^. By combining and 
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Pi, we found the stellar radius to be 

i?* = 0.210 ± O.OO7i?0 (method A). 



(9) 



Finally, assuming that our transit depths were not af- 
fected by starspots (i.e., Cj = 0) we used the radius ra- 
tio derived from our observations, Rp/Ri, = 0.11610 ± 
0.00048, to compute the planetary radius: 



i?p = 2.650 ± 0.090 i?e (method A). 



(10) 



This is in agreement with that found by Charbonneau et 
al. (2009), who found Rp = 2.678±0.130 i?©. Taking the 
more conservative stance that spots may have biased the 
radius ratio measurement, we may set an upper bound 
on the radius of the planet: Rp < 2.81 Rq with 95% 
confidence. 

6.2. Method B: Stellar evolutionary models 

The mass-luminosity relationship of Delfosse et 
al. (2000) has the virtue of being highly empirical, as it 
is based on observations of astrometrically resolved, de- 
tached eclipsing binaries with measured parallaxes. How- 
ever, with only 16 systems as inputs, the empirical rela- 
tion cannot be expected to account for all the relevant 
variables such as age and metallicity. An alternate ap- 
proach is to trust theoretical models of stellar evolution, 
which predict the radius of a star as a function of its 
mass, age, and metallicity, and are calibrated by obser- 
vations when possible. 

To try this approach, we used the stellar evolution- 
ary models of Baraffe et al. (1998), which are presented 
as a series of stellar parameters as a function of age 
(isochrones) for a given metallicity. Since the metallicity 
of GJ 1214 is unknown we assumed a solar metallicity. 
We assigned a likelihood £(j) oc exp[— x^(j)/2] to every 
isochrone entry j such that 



E 



(11) 



where {pk} is the set of stellar parameters that is sub- 
ject to observational constraints (e.g. absolute magni- 
tudes in certain bandpasses), ap^, are the corresponding 
1(7 uncertainties, and Apk arc the differences between the 
isochrone values and the "observed" values. Then, the 
most likely value for a given stellar parameter p (such as 
stellar radius) was found by taking a weighted average 
over all the isochrone points. 



(12) 



The uncertainty in this parameter is taken to be 
the square root of the weighted variance, i.e., CTp = 

y/ (p^) — where averages are taken as defined aboveQ 
We used this technique with constraints based on the 
observed apparent magnitudes in the J, H, and K 
bandpasses and the parallax, all from Charbonneau et 

^ This approach is closely related to the technique advocated by 
Torres et al. (2008). The relatively coarse sampling of the Baraffe 
et al. (1998) isochrones precluded the use of more sophisticated 
numerical integration techniques such as those described by Carter 
et al. (2009). 



al. (2009). The results were A/* 0.131 ±0.044 Mq and 
Ri, =0.191 ±0.019 i?0. 

This method gives a significantly lower mass than the 
first method. This difference is almost entirely due to 
the additional degree of freedom of the stellar age: the 
theoretical evolutionary tracks allow for GJ 1214 to have 
any age, while the Delfosse et al. (2000) relation is a 
consensus result based on a fit to stars of varying age, 
most of which are probably older than a few Gyr. 

This naturally raises the question: what is known 
about the age of GJ 1214? As noted by Charbonneau 
et al. (2009), and confirmed in this work, the optical 
light variations of GJ 1214 are only a few percent in am- 
plitude, suggestive of a more mature star. Moreover, age 
estimates based on the star's space velocity and the es- 
timated stellar rotation period also favor an older star, 
leading Charbonneau et al. (2009) to estimate an age 
of 3-10 Gyr. To account for this evidence, we repeated 
our isochrone analysis but with a prior constraint on the 
stellar age. Specifically, we used a Gaussian prior on 
logiQ(age) with a median value of 10 and a standard de- 
viation of 0.4. With this extra constraint, the stellar 
mass was found to be M^, = 0.156 ± 0.006 Mq, in excel- 
lent agreement with the results of the first method. 

This substantiates our claim that the main difference 
between the two methods of estimating the stellar mass 
is the treatment of stellar age. However, there remains a 
significant discrepancy in the stellar radius. The evolu- 
tionary models with the "old age" prior give 



Ri, = 0.179 ± O.OO6i?0 (method B). 



(13) 



which is about 15% smaller than the result of the first 
method. Again taking the radius ratio to be Rp/R^ = 
0.11610 ± 0.00048, the planetary radius is 



Rp = 2.270 ± 0.080 i?® (method B). 
6.3. Possible resolutions of the two methods 



(14) 



To summarize the preceding results, we have inves- 
tigated two methods of estimating the stellar mass and 
radius based on the available information. The two meth- 
ods can be made to agree on the stellar mass, but they 
refuse to agree on the stellar radius. The first method is 
based on an empirical mass-luminosity relation, and the 
stellar mean density derived from the transit light curve. 
The second method is based on stellar evolutionary mod- 
els constrained by the infrared absolute magnitudes and 
an assumed age >1 Gyr. The first method gives a stellar 
radius that is larger by 15%, which is approximately 4.5 
times larger than the internal uncertainty in either esti- 
mate. The discrepancy is even starker if one compares 
the mean stellar density derived from the transit light 
curve (24.1 ± 1.7 g cm^"^) with that derived from the 
stellar evolutionary models (38.4 ± 2.1 g cm^'^), as they 
disagree by 7a. 

It has long been appreciated that stellar evolutionary 
models have trouble predicting the masses and radii of 
M dwarfs in detached eclipsing binary systems (see, e.g., 
Hoxie 1973, Lacy 1977), in the sense that the models tend 
to underpredict the observed radius for a given mass by 
about 10% (Ribas 2006). It was for this reason that 
Charbonneau et al. (2009) discounted the results based 
on the evolutionary models, and why subsequent authors 
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have done the same in their discussions of GJ 1214. How- 
ever, it has recently been argued that the faihngs of the 
evolutionary models are confined to the mass range 0.3- 
0.7 Mq, and even more specifically to stars that have 
been tidally spun up (and made more active) by a close 
stellar companion (see, e.g., Torres et al. 2006, Lopez- 
Morales 2007, Morales et al. 2008). Stars below a thresh- 
old mass of «0.32 Mq are expected to be fully convective, 
and seem to be well described by the evolutionary mod- 
els (see, e.g., Lopez-Morales 2007, Demory et al. 2009, 
Vida et al. 2009). To the extent this is true, we would 
not expect GJ 1214 — a single, low-activity star of mass 
<0.3 Mq — to be affected by the problems that plague the 
theories of earlier-type M dwarfs in close binary systems. 

We are thereby motivated to seek alternative resolu- 
tions to the discrepancy between the two methods of es- 
timating the mass and radius of G J 1214. Of course there 
is always the possibility that a key input such as the par- 
allax or infrared magnitude is faulty, but in the sections 
to follow we discuss some possible resolutions in which 
all of the data are taken at face value. 

6.3.1. A metal-poor star? 

We have worked exclusively with solar metallicity 
isochrones. Lower-metallicity isochrones generally pre- 
dict a larger radius for a given mass, in the relevant re- 
gion of parameter space. Schlaufman (2010) have pre- 
sented a simple method to estimate the metallicity of an 
M dwarf based on its observed V and K magnitudes. 
For GJ 1214, his method gives [Fe/H]= —0.16, suggest- 
ing that the star is only moderately metal-poor. This 
value for the metallicity would not affect the theoreti- 
cally predicted radius by enough to resolve the discrep- 
ancy. Lopez-Morales (2007) found that the Baraffc et 
al. (1998) isochrones predict a variation in stellar ra- 
dius of only about 3% for metallicities ranging from 
0.0 to —0.5. Moreover, for a sample of low-mass stars 
with -0.5 <[Fe/H]< 0, Demory et al. (2009) showed 
that there is no significant correlation between measured 
metallicity and the fractional difference between the mea- 
sured stellar radius and that found assuming solar metal- 
licity. Therefore, while it is possible that GJ 1214 is a 
metal-poor star, we do not consider this to offer a likely 
resolution of the discrepancy we have noted between the 
two methods of estimating the stellar radius. 

6.3.2. A young star? 

If the star were relatively young and still contracting 
onto the main sequence, then the evolutionary models 
would predict a larger radius, relieving the discrepancy. 
To investigate this possibility we repeated our isochrone 
analysis, but this time with a flat prior on the stellar 
age and a Gaussian prior on the stellar mean density to 
conform with the transit light curve analysis. The result 
was that the age of the star must be about 100 Myr. 
This would conflict with the evidence for an older age, 
namely, the observed lack of strong chromospheric activ- 
ity and the kinematics (see § 16. 2p . There is also the low 
probability that we would happen to observe this star at 
such an early phase of its long life. For these reasons, a 
young age for GJ 1214 does not seem to be a promising 
solution. 

6.3.3. A star spot- corrupted estimate ofa/R^,? 



It is conceivable that unidentified starspot anomalies 
have significantly biased our estimate of a/i?* from the 
transit light curves (see ? I3.1[) . Could this be responsible 
for the 7cr discrepancy in p-^ between the two methods of 
characterizing the star? If this were the case, then the 
more trustworthy estimate of p^, would be the value from 
the isochrone analysis. 

This seems unlikely, partly because the residuals to our 
transit light curves do not display any anomalies beyond 
the two that we have already identified, and partly be- 
cause in this scenario it would be difficult to understand 
the collection of transit depths. Specifically, we can use 
the value of from the evolutionary models as an in- 
put to Equation ([5|), to derive the planet-to-star radius 
ratio. The resuh is i?p/i?* = 0.15 ± 0.01. This conflicts 
with the upper limit Rp/R^, < 0.1161 that we derived in 
§ O under the assumption that the transit depths are af- 
fected by cool spots on the stellar disk. Thus, one would 
have to suppose that nearly all of the transit depths were 
biased to smaller values by numerous starspot crossings 
throughout the transits. We do not consider a conspir- 
acy of starpot anomalies to be a satisfactory solution to 
the problem. 

6.3.4. Is the orbit eccentric? 

If the orbit of GJ 1214b is not circular, then the param- 
eter that is being determined by the light curve analysis 
is not a/Ri, and our application of Equation ([8]) is er- 
roneous. The correct procedure must take into account 
the speed of the planet at inferior conjunction, which is 
a function of the eccentricity and argument of periccn- 
ter (see, e.g., Eqns. 16 and 27 of Winn 2010 or Kipping 
2010). The end result is that the mean density of the 
planet would be modified as follows 



= P*, circ — : , (.l&j 

\ 1 + e sm Lo I 

where p^,^ circ is the mean density that is calculated under 
the assumption of a circular orbit. Therefore, it might 
be possible to reconcile the two different methods for 
estimating the stellar density, for suitable choices of e 
and w. 

The orbital eccentricity has been assumed to be zero, 
because of the expectation that tidal dissipation has 
damped out any initial eccentricity to 10""^ or below. 
However, the published RV data only permit a coarse up- 
per limit of e<0.27 with 95% confidence (Charbonneau 
et al. 2009). To investigate the possibility of an eccentric 
orbit, we assumed that the isochrone-based estimate of 
the mean stellar density (38.4 ±2.1 g cm^"^) is accurate. 
We then derived the constraints on the orbital eccen- 
tricity that would be required for consistency with the 
transit light curve analysis. The minimum eccentricity is 
obtained for the case ui = 90°, corresponding to transits 
occurring at pericenter. In that case, e = 0.138 ± 0.034. 
Solutions with e < 0.27 can be found for values of uj 
between 25° and 155°. 

The only objection is that tidal dissipation should have 
damped out this eccentricity. The characteristic circu- 
larization timescale is Tc ^ 10 Myr assuming a tidal 
quality factor Qp = 100 for the planet, a rough order- 
of-magnitude estimate for a solid planet. This is much 
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shorter than the estimated stellar age. For an icy body 
like Neptune with Qp ^ 10"* — 10^, the circularization 
timescale would be 1-10 Gyr, which is at least not over- 
whelmingly shorter than the system age. The mecha- 
nisms and timescales for tidal dissipation are not under- 
stood from first principles, and are poorly constrained 
by observations. And it should be kept in mind that all 
of the known transiting Neptune-like planets (GJ 436, 
Kepler-4, and HAT-P-11) all have significantly eccentric 
orbits. For these reasons, of all the possibilities we have 
discussed, we find an eccentric orbit to be the most at- 
tractive solution to the radius discrepancy problem. 

7. ANALYSIS OF MIDTRANSIT TIMES: 
LIMITS ON PERTURBER MASS 

We recalculated the transit ephemeris utilizing 14 of 
the 16 midtransit times given in Table [5] We excluded 
the two midtransit times inferred from partial transit 
light curve data (epochs 195 and 200), out of concern 
that the second-order airmass correction could not be 
determined as well in those cases. We fitted a linear 
function of transit epoch E, 

T,iE) = T,{0) + EP. (16) 

The results were Tc(0) = 2454980.7487955 ± 0.000045 
(BJDtdb) and P = 1.58040482 ± 0.00000024 days and 
the fit had = 19-2 with 12 degrees of freedom. For- 
mally the fit is inconsistent with the linear model with 
90% confidence, but we do not consider this to be a sig- 
nificant detection of timing anomalies. The right panel of 
Fig. [5] shows the residuals to the linear fit. This plot also 
shows four data points from Charbonneau et al. (2009), 
which are consistent with the linear model. 

Since there is no clear evidence for transit timing vari- 
ations (TTVs), we may use the timing data to place up- 
per limits on the mass of a hypothetical second planet 
that would perturb the orbit of the transiting planet. 
To do so, we used an implementation of the algorithm 
advocated by Steffen & Agol (2005). We explored the 
parameter space of the hypothetical perturber's mass 
and orbital period and phase. We assumed the orbits 
of the perturber and GJ1214 were circular and coplanar. 
Regions in parameter space yielding dynamical instabil- 
ity, determined following the prescription by Barnes & 
Greenberg (2006), were not sampled. For each orbital 
period ranging from 0.3-8 days and sampling all orbital 
phases, the mass of the perturber was increased until 
the computed transit times would fit our data Ax^ = 9 
worse than a linear ephemeris; this mass is interpreted 
as a 99.7% confidence upper-limit. 

Fig. [8] shows the constraints on the perturber mass 
as a function of period ratio, as determined from this 
analysis. For reference, on the right-hand axis we have 
included the masses of GJ 1214b (6.6 M^), the Moon 
(0.01 A/®), and Eris (0.003 M®). In this plot we have 
also indicated the zones of dynamical instability, and of 
potential habitability (equilibrium temperature between 
273-373 K, for an assumed albedo of zero). We have 
also plotted a line corresponding to an RV amplitude of 
2 m s~^, probably the best achievable detection limit 
based only on RV observations for the near future. 

The mass constraints on the perturber are more re- 
strictive near the mean-motion resonances and most re- 
strictive at the low-order resonances, particularly for the 



interior and exterior 2:1 resonances. For example, a per- 
turber at the interior 2:1 resonance having mass near 
that of Eris would have induced detectable TTVs with 
the present data. 

8. DISCUSSION 

One of our goals in this study was to improve on the es- 
timates of the basic parameters of GJ 1214b. Yet despite 
having undertaken many high-precision observations of 
transits, we have not significantly improved on the esti- 
mate of the planetary radius. In fact, our analysis has 
led us to conclude that the radius is even more uncertain 
than had been recognized previously. This is for two 
reasons. First, the clear evidence for starspots in our 
most precise light curves has caused us to consider the 
possible effects of stellar activity on the analysis of tran- 
sit photometry. Second, and more importantly, we have 
found a significant disagreement between two methods 
of estimating the stellar and planetary dimensions, with 
no good reason why either one should be disregarded or 
considered less trustworthy. Both of these complications 
were known prior to this work (e.g., Charbonneau et al. 
2009), but we have brought them into sharper focus. 

The problem of star spots can be mitigated by observ- 
ing in longer- wavelength bandpasses. This is because the 
flux contrast between two blackbodies of different tem- 
peratures is a decreasing function of wavelength. For 
example, if we assume the particular starspots that pro- 
duced anomalies in our data are approximately the same 
size as the planet, then they are «150 K cooler than 
the stellar photosphere (see 13.21 for details). This cor- 
responds to a surface-brightness ratio of about 0.67 be- 
tween the spots and the surrounding photosphere, at an 
observing wavelength of 0.6 /im. At 3.5 /im, the surface- 
brightness ratio would be about 0.91, representing a 
smaller contrast and correspondingly smaller starspot- 
induced effects. 

The problem of the stellar radius will be more difficult 
to solve, and there is much at stake. The mean planetary 
density could be 1.89±0.33 g cm'^ or 3.03±0.50 g cm'^, 
depending on which route is taken to estimate the stellar 
radius. The lower density would imply that the planet 
must have a dense gaseous atmosphere, for which there 
are many intriguing possible origins and compositions 
(Rogers & Seager 2010, Mfller-Ricci & Fortney 2010). 
In contrast, the higher density could be consistent with 
a solid planet with a very thin (or nonexistent) atmo- 
sphere. We have discussed several possible resolutions 
of this discrepancy, and argued that the most attrac- 
tive possibility is that the planet has an eccentric orbit, 
e w 0.14. This hypothesis can be tested by gathering 
additional RV data, or by measuring the time of occul- 
tation of the planet by the star. Neither task will be 
easy. Assuming the minimum detectable eccentricity to 
vary inversely as the square root of the number of RV 
data points, one would need approximately 4 times as 
many data points (with the same precision as the exist- 
ing data) for a 2cr detection of an eccentricity of 0.14. 
Likewise, the occultation depth is expected to be smaller 
than 0.1% at 3.5 fiui, the region accessible to the Spitzer 
space telescope, which probably offers the best prospects 
for such a detection. 

This experience with GJ 1214 invites some broader re- 
marks about the suitability of M dwarf stars as targets 
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Fig. 8. — Upper mass limits for a hypothetical pcrturbcr as a function of perturber period normalized to the period of GJ 1214b, Pc/ Pb- 



in surveys for transiting planets. The advantages of M 
dwarfs have been discussed by Nutzman & Charbonneau 
(2009), among others. Their smaher radii and masses al- 
low for smaller planets to be detected, at a given signal- 
to-noise ratio. They are numerous in our Galactic neigh- 
borhood. Transits of planets in the habitable zone in 
particular are more probable and, hence, more frequent 
than they are around higher-mass stars. These are deci- 
sive advantages, fully justifying the ongoing efforts that 
concentrate exclusively on M dwarfs. However, there are 



two disadvantages of low-mass stars. They tend to have 
larger spots and an overall higher level of stellar activity 
than more Sun-like stars, which will interfere with the 
measurement of the basic transit parameters (c.f. Seager 
& Deming 2009). And, it will be harder to obtain re- 
liable estimates of the stellar mass and radius, because 
of possible limitations of stellar evolutionary models in 
the range 0.3-0.7 M© and, more generally, because of 
our more limited empirical knowledge of the lowest-mass 
stars. 
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